Generalized Additive Models
Teodoro MOUNIER TEBAS & Diane THIERRY
11/25/2020
Sommaire
- Introduction
- Manipulations et traitement des données
- Modélisations des temps des jobs
- Soumission des prévisions et conclusion
Introduction
Le leader mondial ‘Talend’ de l’intégration et de l’intégrité des données dans le cloud, a annoncé le 30 juillet 2019 la disponibilité de Stitch Data Loader sur la marketplace d’Amazon Web Services (AWS). Déjà pertinente par sa portée et rendue d’autant plus accessible depuis cette annonce, la plateforme Stitch est devenue un incontournable dans la procuration de données ouvertes. Stitch Data Loader facilite en effet l’ingestion de données pour l’analyse et s’adapte automatiquement à l’augmentation du volume de données.
Cette plateforme est une composante de Talend Data Fabric, qui est une suite d’applications cloud conçues pour aider les clients à collecter, gérer, transformer et partager données et informations.
Pour ce projet, ont été récupérées différentes données issues de la plateforme Stitch concernant différentes requêtes effectuées, appelées « jobs » dont nous chercherons à prédire la durée en secondes.
- Nous disposons pour cela de 2 bases :
- une base d’entraînement (« train ») composée de 9.100 observations et 5 variables
- une base de test (« test ») composée de 1.700 observations et 5 variables
Le but de ce projet est de créer des modèles additifs généralisés qui prédisent au mieux le temps de durée de chaque job sur la plateforme Stitch. Nous construirons pour cela différents GAM sur la base d’entraînement, modèles que nous chercherons à optimiser pour les rendre les plus performants possible, puis nous les appliquerons sur la base ‘test’ où nous espérons les erreurs les plus faibles.
Dans un premier temps nous appelons les librairies qu’il nous faut pour réaliser le projet :
library(tidyverse)
library(readr)
library(gganimate)
library(viridis)
library(gridExtra)
library(summarytools)
library(kableExtra)
library(knitr)
library(flextable)
library(tsoutliers)
library(EnvStats)
library(stringr)
library(seastests)
#install.packages("rAmCharts")
library(rAmCharts)
library(mgcv)
library(car)1-Manipulations et traitement des données
1.1 Téléchargement des données
Nous avons téléchargé les données depuis la plateforme Kaggle :
train <- read_csv('./Données/data_train.csv')
test <- read_csv('./Données/data_test.csv')Comme évoqué précédemment, les données brutes contiennent 5 variables et 9100 observations.
Nous disposons initialement de 6 variables :
- JOB_DURATION : durée de la requête Stitch, en secondes - notre variable à expliquer (\(Y_i\)) qui n’est disponible qu’en base ‘train’ puisque nous chercherons à la prédire avec la base test
- TAP_EXIT_STATUS : variable binaire qui indique si le job s’est exécuté, 0 si oui et 1 sinon
- START_TIME : l’heure de début de la requête
- CONNECTION_ID : variable composée de 3 modalités indiquant l’identifiant de la chronique (1, 2 ou 3)
- TAP_VERSION : version du Terminal Access Point utilisée qui est plus ou moins récente, le TAP est un dispositif permettant de surveiller un réseau informatique sans le perturber.
- ID : disponible uniquement dans la base ‘test’, permet d’identifier l’utilisateur
1.2 Premières informations sur les variables
Nous regardons les 5 premières lignes de notre jeu de données pour nous rendre compte des informations qu’il contient. Nous pouvons aussi regarder le format des variables grâce à la commande str(). Puis, le summary de la base nous donnera un aperçu des distributions des variables avec les valeurs minimales et maximales de chacune d’elle ainsi que la moyenne, les quartiles etc. Le but de cette analyse exploratoire étant de cerner au mieux les données pour éventuellement créer de nouvelles variables ou transformer les existantes.
| CONNECTION_ID | START_TIME | TAP_VERSION | TAP_EXIT_STATUS | JOB_DURATION |
|---|---|---|---|---|
| 1 | 2019-10-11 13:41:00 | 1.15.1 | 0 | 9 |
| 1 | 2019-10-11 14:41:00 | 1.15.1 | 0 | 9 |
| 1 | 2019-10-11 15:41:00 | 1.15.1 | 0 | 12 |
| 1 | 2019-10-11 16:41:00 | 1.15.1 | 0 | 10 |
| 1 | 2019-10-11 17:41:00 | 1.15.1 | 0 | 11 |
| 1 | 2019-10-11 18:41:00 | 1.15.1 | 0 | 10 |
| 1 | 2019-10-11 19:41:00 | 1.15.1 | 0 | 10 |
| 1 | 2019-10-11 20:41:00 | 1.15.1 | 0 | 11 |
| 1 | 2019-10-11 21:41:00 | 1.15.1 | 0 | 10 |
| 1 | 2019-10-11 22:41:00 | 1.15.1 | 0 | 11 |
On constate d’emblée que notre variable à expliquer “JOB_DURATION” correspond à la cinquième colonne, et que la variable “START_TIME” regroupe à la fois les dates et les heures du début des jobs. À partir des 10 premières observations, on voit que les variables “CONNECTION_ID”, “TAP_VERSION” et “TAP_EXIT_STATUS” restent inchangées, ce qui n’est pas le cas de “JOB_DURATION” qui évolue en fonction du temps récupéré par la variable “START_TIME”. De plus, les durées des 10 premiers jobs sont très courtes, allant de 9 à 12 secondes.
str(train)
## Classes 'spec_tbl_df', 'tbl_df', 'tbl' and 'data.frame': 9100 obs. of 5 variables:
## $ CONNECTION_ID : num 1 1 1 1 1 1 1 1 1 1 ...
## $ START_TIME : POSIXct, format: "2019-10-11 13:41:00" "2019-10-11 14:41:00" ...
## $ TAP_VERSION : chr "1.15.1" "1.15.1" "1.15.1" "1.15.1" ...
## $ TAP_EXIT_STATUS: num 0 0 0 0 0 0 0 0 0 0 ...
## $ JOB_DURATION : num 9 9 12 10 11 10 10 11 10 11 ...
## - attr(*, "spec")=
## .. cols(
## .. CONNECTION_ID = col_double(),
## .. START_TIME = col_datetime(format = ""),
## .. TAP_VERSION = col_character(),
## .. TAP_EXIT_STATUS = col_double(),
## .. JOB_DURATION = col_double()
## .. )Sur nos 5 variables seules 3 sont numériques ; “TAP_VERSION” est en format ‘character’ et “START_TIME” au format ‘POSIXct’ ce que nous devrons donc changer pour un meilleur traitement des données.
summary(train)
## CONNECTION_ID START_TIME TAP_VERSION
## Min. :1.000 Min. :2019-10-11 13:41:00 Length:9100
## 1st Qu.:1.000 1st Qu.:2020-01-11 05:31:00 Class :character
## Median :1.000 Median :2020-02-27 21:50:30 Mode :character
## Mean :1.626 Mean :2020-03-19 21:41:01
## 3rd Qu.:2.000 3rd Qu.:2020-07-13 01:46:30
## Max. :3.000 Max. :2020-08-29 11:09:00
##
## TAP_EXIT_STATUS JOB_DURATION
## Min. :0.00000 Min. : 2.0
## 1st Qu.:0.00000 1st Qu.: 10.0
## Median :0.00000 Median : 14.0
## Mean :0.00099 Mean : 400.5
## 3rd Qu.:0.00000 3rd Qu.: 181.0
## Max. :1.00000 Max. :4349.0
## NA's :11On voit que sur les 3 modalités de connexion, la 1 regroupe au moins la moitié des observations puisque la médiane est de 1. On constate aussi que les requêtes de la base train s’étendent du 11 octobre 2019 au 29 août 2020. Concernant la variable “TAP_EXIT_STATUS” c’est la seule qui contient des valeurs manquantes, il faudra donc retirer ces 11 observations pour les modélisations. Enfin, on voit que les durées des requêtes vont de 2 secondes à 4349 soit 1 heure 12 minutes et 29 secondes, ce qui est considérable. On peut donc supposer que le job n’a pas abouti (TAP_EXIT_STATUS=1) ou bien qu’il existe une grande disparité des temps au sein des périodes. On voit ainsi que la moyenne est de 400 secondes tandis que la médiane est de seulement 14.
1.3 Visulation des données
Avant toute transformation des variables nous allons visualiser la variable à expliquer en fonction du temps.
# On plot la variable à expliquer fonction du temps
plot(train$JOB_DURATION ~ train$START_TIME, type='l', xlab="Temps", ylab="Durées des jobs en secondes", main="Évolution de la durée des requêtes dans le temps")
Au vu du graphique ci-dessus on peut constater qu’il y a 3 périodes correspondant aux 3 chroniques :
- la première s’étend d’octobre à mai 2020 avec une volatilité très faible donc un temps des jobs très court
- la seconde de juillet à septembre comprenant donc l’été, avec des temps légèrement plus longs que ceux de la première période liés probablement aux vacances scolaires d’été
- la dernière s’étend de janvier à fin mars avec cette fois des temps beaucoup plus importants qui semblent varier autour de 2000 secondes soit 33 minutes
La connexion ID n°3 se distingue donc largement des chroniques 1 et 2 par sa grande volatilité de temps des jobs, malgré le fait qu’elle ait été capturée sur une sous-période de la chronique 1. Compte tenu du fait que les valeurs des jobs diffèrent beaucoup d’une chronique à l’autre, nous allons procéder à la division de notre base en 3 sous bases pour chaque chronique, dans le but d’observer plus précisément les variations des temps des jobs.
ID 1
# On créé une sous-base pour la chronique n°1
df_ID1 <- train %>% filter(CONNECTION_ID==1)
summary(df_ID1)
## CONNECTION_ID START_TIME TAP_VERSION TAP_EXIT_STATUS
## Min. :1 Min. :2019-10-11 13:41:00 Length:5000 Min. :0
## 1st Qu.:1 1st Qu.:2019-12-02 14:26:00 Class :character 1st Qu.:0
## Median :1 Median :2020-01-23 16:11:00 Mode :character Median :0
## Mean :1 Mean :2020-01-23 22:53:32 Mean :0
## 3rd Qu.:1 3rd Qu.:2020-03-16 09:56:00 3rd Qu.:0
## Max. :1 Max. :2020-05-07 12:41:00 Max. :0
## JOB_DURATION
## Min. : 7.00
## 1st Qu.: 10.00
## Median : 10.00
## Mean : 12.25
## 3rd Qu.: 12.00
## Max. :128.00La sous base de la chronique 1 ainsi créée contient 55% de la base train pour un total de 5000 observations. Nous constatons que pour la variable “TAP_EXIT_STATUS” il n’y a plus de valeur manquante, celles-ci devaient donc concerner les chroniques 2 ou 3 (nous verrons cela juste après). En revanche si l’on regarde la durée des jobs, on voit que le temps maximal est de 128 secondes, il s’agit sûrement d’une valeur atypique qui peut gêner nos modélisations et que nous allons donc traiter dès maintenant.
amBoxplot(df_ID1$JOB_DURATION, xlab="CONNECTION ID n°1", ylab="Durée des jobs", main="Boxplot 'JOB_DURATION' pour la chronique n°1")À partir du graphique ci-dessus, nous observons d’emblée que de nombreuses observations sont supérieures au troisième quantile et peuvent donc être considérées comme atypiques. Cependant, étant donné le nombre important de points concernés il peut s’agir d’un comportement explicable et qui peut apporter une information importante pour les futures modélisations. Nous considérons donc comme atypiques les valeurs qui dépassent 1 minute car il n’y en a que 8 qui sont assez espacées, et procédons au test de Rosner pour confirmer leur atypicité.
rosnerTest(df_ID1$JOB_DURATION, k = 8, alpha = 0.05) # 8 outliers : 128,126,121,114,108,76,75,65# On regarde à quelles dates les points correspondent pour voir si pattern
a <- df_ID1 %>% slice(1423,846,967,2479,487,1711,366,2489)
kbl(a, align = "c") %>% kable_classic(full_width = F, html_font = "Cambria")| CONNECTION_ID | START_TIME | TAP_VERSION | TAP_EXIT_STATUS | JOB_DURATION |
|---|---|---|---|---|
| 1 | 2019-12-09 18:41:00 | 1.15.1 | 0 | 128 |
| 1 | 2019-11-15 17:41:00 | 1.15.1 | 0 | 126 |
| 1 | 2019-11-20 18:41:00 | 1.15.1 | 0 | 121 |
| 1 | 2020-01-22 18:41:00 | 1.15.1 | 0 | 114 |
| 1 | 2019-10-31 18:41:00 | 1.15.1 | 0 | 108 |
| 1 | 2019-12-21 18:41:00 | 1.15.1 | 0 | 76 |
| 1 | 2019-10-26 18:41:00 | 1.15.1 | 0 | 75 |
| 1 | 2020-01-23 04:41:00 | 1.15.1 | 0 | 65 |
D’après le test de Rosner on voit que les 8 observations étaient effectivement atypiques, nous trouvions intéressant de regarder de quelles observations il s’agissait pour noter d’éventuels comportements communs par exemple par heure ou par jour de la semaine. Nous observons alors d’après la table que la plupart des observations atypiques se situent en fin d’après-midi à 18h41 (c’est le cas de 6 observations sur 8). On ne note cependant pas de régularité dans les jours de la semaine puisque l’on retrouve 2 fois mercredi, jeudi et samedi, et 1 fois lundi et vendredi.
# on crée une nouvelle base sans ces points
df_ID1=df_ID1[-c(1423,846,967,2479,487,1711,366,2489),]
dim(df_ID1)
## [1] 4992 5
amBoxplot(df_ID1$JOB_DURATION, xlab="CONNECTION ID n°1", ylab="Durée des jobs", main="Boxplot 'JOB_DURATION' débarrassé des points atypiques") # Graphiques CONNECTION_ID 1
par(mfrow=c(2,2))
plot(JOB_DURATION ~ START_TIME, type='l', data=df_ID1, main="Toute la période", xlab="Temps", ylab="Durées des jobs") # le tout
plot(JOB_DURATION ~ START_TIME, type='l', data=df_ID1[84:107,], main="Un mardi", xlab="Temps", ylab="Durées des jobs") # mardi
plot(JOB_DURATION ~ START_TIME, type='l', data=df_ID1[564:731,], main="Une semaine", xlab="Temps", ylab="Durées des jobs") # une semaine
plot(JOB_DURATION ~ START_TIME, type='l', data=df_ID1[491:1208,], main="Un mois (novembre)", xlab="Temps", ylab="Durées des jobs") # un mois
On peut constater à partir des 4 graphiques des tendances :
- chaque journée commence avec un pic significatif autour de 6 heures du matin, puis un autre plus faible en fin de journée
- ces pics où les temps des jobs sont élevés, sont visibles pour chaque jour de la semaine peu importe qu’il s’agisse d’un week end ou d’un début de semaine
- pour le mois de novembre on retrouve ces pics avec quelques valeurs plus grandes mais non régulières
Évolution de Y sur plusieurs semaines
par(mfrow=c(2,2))
plot(JOB_DURATION ~ START_TIME, type='l', data=df_ID1[60:227,], main="Une semaine", xlab="Temps", ylab="Durées des jobs")
plot(JOB_DURATION ~ START_TIME, type='l', data=df_ID1[228:395,], main="Une semaine", xlab="Temps", ylab="Durées des jobs")
plot(JOB_DURATION ~ START_TIME, type='l', data=df_ID1[396:563,], main="Une semaine", xlab="Temps", ylab="Durées des jobs")
plot(JOB_DURATION ~ START_TIME, type='l', data=df_ID1[564:731,], main="Une semaine", xlab="Temps", ylab="Durées des jobs") 
Évolution de Y sur plusieurs mois
par(mfrow=c(2,2))
plot(JOB_DURATION ~ START_TIME, type='l', data=df_ID1[491:1208,], main="Un mois (novembre)", xlab="Temps", ylab="Durées des jobs")
plot(JOB_DURATION ~ START_TIME, type='l', data=df_ID1[1209:1926,], main="Un mois (décembre)", xlab="Temps", ylab="Durées des jobs")
plot(JOB_DURATION ~ START_TIME, type='l', data=df_ID1[1927:2644,], main="Un mois (janvier)", xlab="Temps", ylab="Durées des jobs")
plot(JOB_DURATION ~ START_TIME, type='l', data=df_ID1[2645:3362,], main="Un mois (février)", xlab="Temps", ylab="Durées des jobs") 
ID 2
# On créé une sous-base pour la chronique n°2
df_ID2 <- train %>% filter(CONNECTION_ID==2)
summary(df_ID2)
## CONNECTION_ID START_TIME TAP_VERSION
## Min. :2 Min. :2020-07-08 09:39:00 Length:2500
## 1st Qu.:2 1st Qu.:2020-07-21 10:01:30 Class :character
## Median :2 Median :2020-08-03 10:24:00 Mode :character
## Mean :2 Mean :2020-08-03 10:24:00
## 3rd Qu.:2 3rd Qu.:2020-08-16 10:46:30
## Max. :2 Max. :2020-08-29 11:09:00
##
## TAP_EXIT_STATUS JOB_DURATION
## Min. :0.000000 Min. : 67.0
## 1st Qu.:0.000000 1st Qu.:165.8
## Median :0.000000 Median :172.0
## Mean :0.001601 Mean :178.4
## 3rd Qu.:0.000000 3rd Qu.:182.0
## Max. :1.000000 Max. :688.0
## NA's :1Nous réitérons l’analyse pour la chronique 2 en créant une autre sous base pour cet ID qui contient cette fois 2500 observations. On voit cette fois que la durée moyenne des jobs est beaucoup plus importante que pour la chronique 1 variant autour de 178 secondes. De plus nous observons qu’il y a une valeur manquante pour la variable “TAP_EXIT_STATUS” que nous retirons avant la détection des outliers.
#on supprime la valeur manquante pour "TAP_EXIT_STATUS" (6 aout 2020)
df_ID2 <- na.omit(df_ID2)
#on regarde le boxplot
amBoxplot(df_ID2$JOB_DURATION, xlab="CONNECTION ID n°2", ylab="Durée des jobs", main="Boxplot 'JOB_DURATION' pour la chronique n°2")| CONNECTION_ID | START_TIME | TAP_VERSION | TAP_EXIT_STATUS | JOB_DURATION |
|---|---|---|---|---|
| 2 | 2020-07-30 12:09:00 | 1.4.33 | 0 | 688 |
| 2 | 2020-07-28 17:39:00 | 1.4.33 | 0 | 604 |
Comme visible sur le boxplot il y a 2 valeurs extrêmes correspondant aux observations 1062 et 977 (détectées avec le test de Rosner) ; les 2 durées atypiques d’environ 10 minutes d’attente ont eu lieu à deux journées d’intervalle fin juillet 2019. Les temps de la variable “JOB_DURATION” sont maintenant compris entre 151 et 396 secondes. En moyenne le temps d’attente pour la chronique 2 est 15 fois plus long que pour la chronique 1.
Par rapport aux graphiques de la chronique une, on retrouve le pic du temps des jobs de début de journée assez grand et un pic plus faible en fin d’après-midi. En revanche on constate que le pic de début de journée est légèrement plus tard qu’en période “scolaire” puisque nous le rappelons la chronique 2 concerne les mois de juillet et août 2020.
Évolution de Y sur plusieurs semaines
ID 3
# On créé une sous-base pour la chronique n°3
df_ID3 <- train %>% filter(CONNECTION_ID==3)
summary(df_ID3)
## CONNECTION_ID START_TIME TAP_VERSION
## Min. :3 Min. :2020-01-08 03:02:00 Length:1600
## 1st Qu.:3 1st Qu.:2020-01-24 17:45:00 Class :character
## Median :3 Median :2020-02-10 09:30:00 Mode :character
## Mean :3 Mean :2020-02-10 10:02:12
## 3rd Qu.:3 3rd Qu.:2020-02-27 01:15:00
## Max. :3 Max. :2020-03-14 19:00:00
##
## TAP_EXIT_STATUS JOB_DURATION
## Min. :0.000000 Min. : 2
## 1st Qu.:0.000000 1st Qu.:1664
## Median :0.000000 Median :1824
## Mean :0.003145 Mean :1961
## 3rd Qu.:0.000000 3rd Qu.:2165
## Max. :1.000000 Max. :4349
## NA's :10Finalement, nous réalisons une troisième fois l’analyse pour la dernière connexion : la base contient alors 1600 observations et 10 valeurs manquantes que nous enlevons de la base. Pour cette période de janvier à mars 2020 on voit que les durées des requêtes sont encore plus élevées que la chronique 2.
| CONNECTION_ID | START_TIME | TAP_VERSION | TAP_EXIT_STATUS | JOB_DURATION |
|---|---|---|---|---|
| 3 | 2020-03-05 21:00:00 | 2.6.5 | 0 | 4349 |
| 3 | 2020-01-08 06:01:00 | 2.6.4 | 0 | 660 |
| 3 | 2020-01-16 05:01:00 | 2.6.4 | 1 | 841 |
| 3 | 2020-01-21 20:00:00 | 2.6.4 | 0 | 562 |
| 3 | 2020-02-09 07:01:00 | 2.6.4 | 1 | 16 |
| 3 | 2020-03-02 01:01:00 | 2.6.5 | 1 | 7 |
| 3 | 2020-03-13 01:01:00 | 2.6.5 | 1 | 8 |
À partir du boxplot on voit que les observations qui semblent atypiques sont celles qui ont un temps d’éxécution inférieur à 1000 et supérieur à 3000 secondes, nous décidons pour notre part de retirer les valeurs non comprises dans l’intervalle [1000;4000] où 7 valeurs sont concernées. Il y a ainsi 3 points en janvier, 1 en février et 3 autres en mars - valeurs que nous retirons de la base df_ID3. En moyenne sur cette sous-période de la chronique une, les jobs durent 1969 secondes, soit plus de 30 minutes.
On peut constater plusieurs choses cette fois-ci :
- dans la journée le pic où le temps de job est le plus élevé est en fin d’après midi et non le matin
- dans une semaine on voit la plus grande rapidité des requêtes le week-end hormis le dimanche soir
- les graphiques du mois et de toute la période confirment cet “effet week-end”
Évolution de Y sur plusieurs semaines
Cette partie nous a permis de constater la proximité des temps de jobs pour les chroniques 1 et 2 quoiqu’un peu plus longues pour la chronique 2 du fait peut être de l’été. En revanche les temps d’attente des jobs de la chronique 3 sont eux beaucoup plus importants avec des pics à des horaires différents des autres CONNECTION_ID.
1.4 Transformation des variables
Nous allons procéder à la transformation des variables qui nous ont semblées pertinentes dans l’explication du temps des jobs. Ayant découpé l’échantillon ‘train’ en 3 sous-bases selon la chronique, la seule variable qui reste pertinente pour l’explication des jobs est la variable “START_TIME”. A partir de cette dernière nous allons pouvoir extraire plusieurs informations pour expliquer le temps des jobs. En découpant “START_TIME” en nombre d’heures par jour, des “paterns” ont été mis en évidence quelles que soit la chronique. Lorsque nous avons découpé cette variable par semaine pour voir si le jour de la semaine pouvait influencer le temps des jobs, nous nous sommes rendu compte que seule la chronique 3 semblait être affectée par le samedi et dimanche que nous avons appelé précédemment “effet week-end”. Enfin, lorsque l’on a observé Y sur plusieurs mois nous n’avons pas noté de ressemblances ou de tendances d’un mois à l’autre. Dans tous les cas nous faisons le choix d’extraire de la variable “START_TIME” l’heure de la journée qui nous semble être la variable la plus importante, le jour de la semaine et la semaine de l’année, tels que :
- Week : le numéro de la semaine dans l’année, les valeurs iront donc de 1 à 52
- Day : le jour dans la semaine avec ainsi des valeurs allant de 1 à 7
- Hour : l’heure de la journée
# On crée toutes ces variables en séparant chaque info de 'Start_time'
# Semaine de l'année
df_ID1$Week <- format(df_ID1$START_TIME, format = "%U")
df_ID2$Week <- format(df_ID2$START_TIME, format = "%U")
df_ID3$Week <- format(df_ID3$START_TIME, format = "%U")
# Jours de la semaine
#on récupère le nom de la semaine
df_ID1$Day <- format(df_ID1$START_TIME, format = "%A")
df_ID2$Day <- format(df_ID2$START_TIME, format = "%A")
df_ID3$Day <- format(df_ID3$START_TIME, format = "%A")
#on convertit le nom en chiffre
df_ID1$Day <- str_replace_all(df_ID1$Day, c("Lundi"="1", "Mardi"="2", "Mercredi"="3", "Jeudi"="4", "Vendredi"="5", "Samedi"="6", "Dimanche"="7"))
df_ID2$Day <- str_replace_all(df_ID2$Day, c("Lundi"="1", "Mardi"="2", "Mercredi"="3", "Jeudi"="4", "Vendredi"="5", "Samedi"="6", "Dimanche"="7"))
df_ID3$Day <- str_replace_all(df_ID3$Day, c("Lundi"="1", "Mardi"="2", "Mercredi"="3", "Jeudi"="4", "Vendredi"="5", "Samedi"="6", "Dimanche"="7"))
# Heures
df_ID1$Hour <- format(df_ID1$START_TIME, format = "%H")
df_ID2$Hour <- format(df_ID2$START_TIME, format = "%H.%M") #on ajoute les minutes puisqu'il y a 2 obs par heure pour l'ID2
df_ID3$Hour <- format(df_ID3$START_TIME, format = "%H")
# On garde dans une variable juste la date sans l'heure
df_ID1$Date <- as.Date(df_ID1$START_TIME, format = "%m-%d-%Y")
df_ID2$Date <- as.Date(df_ID2$START_TIME, format = "%m-%d-%Y")
df_ID3$Date <- as.Date(df_ID3$START_TIME, format = "%m-%d-%Y")Les commandes ci-dessus permettent donc de créer les différentes variables. Ce processus d’extraction transforme les variables créées en ‘caractère’, or les modèles GAM ne supportent pas ce format pour les modélisations donc nous les convertissons au format numérique par la fonction as.numeric() seulement pour la variable ‘Hour’ qui prend les heures et les minutes pour l’ID2, sinon nous appliquons la fonction as.integer() pour les autres variables puisque les valeurs sont discrètes.
Nous voilà désormais avec 9 variables par sous-base de chronique. L’intérêt de séparer quelques informations contenues dans la variable ‘START_TIME’ se trouve dans les modélisations GAM : le fait de les avoir passées au format numérique permet d’appliquer des fonctions des variables dans les modèles, et ainsi améliorer les prévisions.
1.5 Statistiques descriptives
De nombreux packages sous R offrent des statistiques intéressantes et concrètes pour aider à comprendre les données, tels que sumarytools, DataExplorer et explore. Nous regardons les principales informations de notre variable à expliquer grâce à la première librairie, en considérant les temps des jobs par chronique car les écarts sont très importants d’une à l’autre comme nous avons pu le montrer précédemment.
| JOB_D ID1 | JOB_D ID2 | JOB_D ID3 | |
|---|---|---|---|
| Mean | 12.11 | 178.12 | 1973.89 |
| Median | 10.00 | 172.00 | 1825.00 |
| Std.Dev | 7.04 | 22.15 | 437.12 |
| Min | 7.00 | 151.00 | 1100.00 |
| Max | 58.00 | 396.00 | 3622.00 |
Sur les bases propres c’est à dire débarrassées des valeurs atypiques, on constate une fois encore l’intérêt de diviser en 3 bases du fait des écarts importants dans les durées des jobs. Pour pallier à de tels écarts, nous aurions pu standardiser les temps pour être sur la même échelle mais cela aurait compliqué les interprétations et non avons donc choisi de garder 3 bases distinctes.
Avec une semaine prise au hasard parmi les observations de la chronique 1, par défaut avec la fonction ggplot() on voit que l’ajustement est mauvais, d’où l’intérêt de réaliser des modèles GAM où l’on pourra spécifier chaque paramètre de la régression.
2-Modélisations des durées des jobs
Nous allons dans cette partie commencer les modélisations des durées des jobs, en construisant par un modèle simple GAM par chronique puis en essayant d’optimiser les paramètres et ainsi améliorer les prévisions. Pour évaluer la performance des modèles nous regarderons différentes choses :
- R2 : qui représente le pourcentage de la variance expliquée par le modèle (sa qualité d’ajustement)
- Significativité des variables : on ne conserve que les modèles où toutes les variables sont significatives en respectant le critère de parcimonie
- Visualisation graphique : sur un même graphique nous confronterons les valeurs prédites par le modèle aux données observées pour une semaine tirée au hasard et non toute la période pour une meilleure lisibilité
Aussi, pour comparer les modèles entre eux nous regarderons les critères suivants :
- GCV : qui donne la moyenne des erreurs des estimations par validation croisée
- AIC = 2 log(L) − 2p avec L le Log de Vraisemblance du modèle et p le nombre de paramètres
- MAPE : qui donne la moyenne des erreurs absolues des prévisions en pourcentage
Toutes ces informations auront pour but de constater la qualité de nos prévisions et ainsi sélectionner le meilleur modèle, c’est à dire celui qui maximise le \(R^2\) et minimise l’AIC, le GCV et le MAPE.
2.1 Fonctionnement d’un modèle GAM
Introduit dès les années 80, le modèle additif généralisé (GAM) permet de traiter des relations non linéaires que peuvent prendre les variables explicatives, en s’affranchissant de la linéarité pour mieux prendre en compte les seuils et les valeurs extrêmes. Le principe d’un GAM est ainsi d’expliquer une phénomène non par des variables explicatives mais par des fonctions de ces dernières : appelées fonctions de lissage non paramétriques ou non spécifiques. Effectivement, aucun paramètre n’est nécessaire à la construction de GAM ; la distribution de la variable à expliquer peut alors prendre la forme d’une loi binomiale, d’une loi de Poisson, de Gamma etc. Le modèle linéaire général est ainsi un cas particulier du modèle linéaire généralisé, où Y suit une loi normale et où les fonctions de liaisons avec les variables explicatives sont de simples fonctions identité telles que f(x)=x.
Pour étudier le comportement des fonctions sur Y, il convient de définir et estimer les fonctions non linéaires \(f_{j}\) qui remplacent les coefficients \(\beta_{j}\). Valables pour les variables quantitatives, ces fonctions permettent de maximiser la qualité de la prévision de Y et ainsi augmenter les performances des modèles construits. Dans le cadre de ce travail nous cherchons à obtenir les meilleures prévisions des durées de différents “jobs” en secondes, c’est pourquoi nous aurons recours aux modèles additifs généralisés.
Plusieurs paramètres entrent en compte lors de l’estimation de GAM et leur ajustement peut permettre une meilleure qualité prévisionnelle, nous jouerons donc sur les paramètres suivants :
- k : le nombre de fonctions de base, plus il augmente plus le nombre d’observations dans la base doit être grand du fait des nombreux paramètres à estimer, autremit dit k découpe l’espace et définit la précision de la régression
- sp : le nombre de splines, une spline représente un polynome de degré \(k\) qui a un instant t vaut 1 et à tous les autres vaut 0 ; on peut alors fixer \(sp\) le nombre maximum de splines sur un espace k. Attention, \(sp=\frac{1}{nb de splines}\) donc plus \(s\) est petit plus le nombre de splines est élevé. Il existe 2 types de splines ; les paramétriques (splines de régression/pénalisés) et les non-paramétriques (splines de lissage)
- method : la méthode d’estimation des paramètres de lissage, lorsque l’on fixe method=“REML” alors le modèle cherche lui-même les meilleurs paramètres de régression. La méthode “GCV” prédit Y par cross-validation et peut aussi être intéressante
- bs : le type de spline de régression parmi “thin plate” (par défaut), “cubic regression”, “cyclic”, “duchon”, “on the sphere” splines etc.
- knots : les nœuds représentent les intervalles auxquels doivent être construites les fonctions de base, on peut donc définir ces intervalles
Par la combinaison de tous ces paramètres on cherche la fonction la plus lisse mais qui soit la plus précise possible pour modéliser Y, on cherche le compromis entre le log de vraisemblance du modèle et sa complexité - compromis définit par le paramètre de lissage \(\lambda\). L’idée d’un GAM est d’utiliser une forme plus flexible de régression entre Y et X ; mais un \(\lambda\) trop petit conduit à une régression trop précise et donc inefficace qu’on appelle over-fiting ou sur-apprentissage (le cas extrême correspondant à une fonction de base par observation). La complexité est aussi observable par l’edf (estimation du nombre de degrés de liberté) qui représente de combien la variable est lissée donc plus l’edf est grand et plus les splines ont une forme complexe.
2.2 Application de modèles GAM
ID 1
Pour les modélisations nous procéderons à la même approche pour les 3 ID ; nous modéliserons dans un premier temps un modèle additif généralisé le plus simple possible, en incluant les variables qui vont avec la chronique mais sans préciser de paramètres énumérés ci-dessus. Puis, petit à petit nous effectuerons une recherche du meilleur modèle en testant différents paramètres jusqu’à trouver les prévisions les plus justes.
Modélisations
1er GAM
# On estime un premier modèle gam le plus simple possible c-à-d. sans rien définir
gam1_ID1 <- gam(formula=JOB_DURATION ~ TAP_VERSION + s(Week) + s(Day,k=1) + s(Hour), data = df_ID1, method="REML")
summary(gam1_ID1) #R2=0.348
##
## Family: gaussian
## Link function: identity
##
## Formula:
## JOB_DURATION ~ TAP_VERSION + s(Week) + s(Day, k = 1) + s(Hour)
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 12.46087 0.15400 80.915 <2e-16 ***
## TAP_VERSION1.16.0 -0.23851 0.46818 -0.509 0.6105
## TAP_VERSION1.17.0 0.07533 1.34894 0.056 0.9555
## TAP_VERSION1.17.1 0.15173 1.07576 0.141 0.8878
## TAP_VERSION1.17.2 -0.17043 0.56229 -0.303 0.7618
## TAP_VERSION1.17.3 -1.64044 0.50608 -3.241 0.0012 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(Week) 3.350 4.064 6.596 2.21e-05 ***
## s(Day) 1.868 1.982 3.506 0.025 *
## s(Hour) 8.975 9.000 290.297 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.348 Deviance explained = 35.1%
## -REML = 15794 Scale est. = 32.355 n = 4992Pour ce premier GAM nous incluons les variables créées ainsi que ‘TAP_VERSION’ de la base de la chronique 1 (‘df_ID1’). On voit alors que cette dernière n’est pas très pertinente puisque seule 1 de ses modalités est significative, au seuil de risque de 5%. Par ailleurs, la variable des heures est elle très pertinente et considérer sa relation à Y de manière non linéaire permet effectivement d’améliorer le modèle (p-value de <2e-16 inférieure à 0.05). Le constat est le même pour les fonctions des variables ‘Week’ et ‘Day’ qui sont elles aussi significatives, cependant on voit que le degré de liberté des jours dans le semaine est très proche de 1 signifiant ainsi une relation presque linéaire avec Y.
plot(gam1_ID1,residuals=FALSE, shade=TRUE, pages=1)
L’utilité de considérer des fonctions des variables explicatives peut se confirmer par les plots des résidus où nous observons clairement des seuils et des variations pour la variable ‘Hour’, une relation légèrement incurvée pour les numéros des semaines, et une quasiment linéaire pour le jour de la semaine sur les temps d’exécution des jobs.
# On plote les valeurs prédites et les valeurs observées pour voir la qualité du modèle simple
fitted <- gam1_ID1$fitted.values
real_ID1 <- df_ID1$JOB_DURATION
date_ID1 <- df_ID1[,2]
prev <- cbind(date_ID1, real_ID1, fitted)
prev <- data.frame(prev)
ggplot(data=prev[60:227,])+
geom_line(mapping=aes(y=real_ID1,x= START_TIME,color="valeurs observées")) +
geom_line(mapping=aes(y=fitted,x= START_TIME,color="valeurs prédites")) +
theme_bw() +
scale_color_manual(values = c(
'valeurs prédites' = 'blue',
'valeurs observées' = 'red')) +
labs(color = 'JOB DURATION', x = "Temps", y = "Durée des jobs",
title = "Valeurs prédites vs. valeurs observées du GAM n°1")
Sur le graphique ci-dessus nous pouvons constater les valeurs prédites par le modèle (en bleu) et les valeurs observées de Y (en rouge) ; nous voyons que ce premier GAM a du mal à modéliser une si grande volatilité et capte mal les valeurs importantes que nous voyons tous les jours à 5h41 pour la chronique 1. On cherche donc a améliorer ce modèle en prenant en compte les constats faits lors de son analyse.
2è GAM
Pour le premier modèle additif généralisé nous avons décidé de ne préciser aucun paramètre pour les fonctions des variables explicatives que nous avons incluses. Dans une deuxième estimation nous définissons le nombre de fonctions de base k par le nombre de valeurs que prennent les différentes variables pour avoir une régression la plus précise possible. Ainsi pour la fonction de l’heure nous posons k=24 puisque les valeurs s’étendent de 0 à 23 heures, de la même manière pour la variable ‘Day’ qui prend des valeurs entre 1 et 7 pour les jours de la semaine nous définissons k=7. Pour la variable ‘Week’ qui indique le numéro de la semaine dans l’année ses paramètres étaient plus difficiles à définir puisque le jeu de données de la chronique 1 est à cheval sur 2 années ; les valeurs de la semaine vont donc de 1 à 52 mais elle ne prend que 30 valeurs différentes (40 à 52è semaine de 2019 et 1 à 18è de 2020). Concernant le type de spline appliqué pour les heures, un comportement cyclique est clairement identifiable donc nous avons ajouté l’argument bs=“cc”, et pour capter la saisonnalité des fonctions de ‘Week’ et ‘Day’ nous avons spécifié bs=“ps” qui indique des p-splines après avoir testé différents types et que celui-ci se soit avéré le plus pertinent.
# Deuxième GAM : on précise les paramètres de celles qui sont utiles
gam2_ID1 <- gam(formula=JOB_DURATION ~ TAP_VERSION + s(Week, k=30, bs="ps") + s(Day, k=7, bs="ps") + s(Hour, k=24, bs="cc"), data = df_ID1, family="gaussian")
summary(gam2_ID1)$r.sq #0.562
## [1] 0.5616597
summary(gam2_ID1)$s.table #2 fonctions significatives
## edf Ref.df F p-value
## s(Week) 3.441082 4.174156 9.900653 5.402969e-08
## s(Day) 2.134445 2.555667 4.123971 9.765975e-03
## s(Hour) 21.928767 22.000000 287.964782 0.000000e+00
plot(gam2_ID1,residuals=FALSE, shade=TRUE, pages=1) # relation quasi linéaire pour 'Day' mais pas pour 'Hour' (cf. edf)
Lorsque nous spécifions les paramètres pour chaque fonction de variable dans le modèle, nous voyons que la qualité du modèle croît beaucoup avec un \(R^2\) plus de 20% supérieur à celui du premier GAM et égal à 56.2%. Par la fonction ‘modele$s.table’ nous voyons la significativité des fonctions des \(X_i\) ; elles le sont toutes au seuil de risque de 5%. Nous constatons aussi que l’edf de la variable des jours de la semaine ‘Day’ est légèrement plus élevé que dans notre premier modèle mais reste proche de la parfaite linéarité (edf=1). Regardons l’évolution des valeurs observées et prédites à nouveau par un graphique temporel.
Comme visible sur le graphique ci-dessus, augmenter le nombre de fonctions de base a permis d’améliorer nos prévisions : nous voyons que les valeurs extrêmes sont beaucoup mieux captées par ce nouveau modèle en revanche les petites fluctuations entre les valeurs importantes sont moins bien modélisées et plus lisses. L’idéal serait donc de trouver le modèle qui capte à la fois les valeurs extrèmes et les variations entre celles-ci.
3è GAM
Dans ces troisièmes estimations nous nous intéressons aux effets non-linéaires bivariés où l’on suppose une corrélation entre certaines variables, ces interactions sont importantes dans les séries temporelles à double saisonnalité. Lorsque l’on inclut dans une même fonction différentes variables explicatives, on suppose que leurs effets se recoupent. Dans un premier temps nous utiliserons une fonction simple d’interaction s(\(X_1\),\(X_2\)) pour voir quelles variables peuvent avoir un effet bivarié, puis on cherchera à préciser cet effet. On ajoute donc au modèle n°2 existant les effets bivariés pour les variables ‘Hour’, ‘Day’ et ‘Week’ en essayant les 3 combinaisons possibles.
# GAM 3 : on ajoute une fonction qui prend un effet bivarié
gam3_ID1 <- gam(JOB_DURATION ~ TAP_VERSION + s(Week, k=30, bs="ps") + s(Day, k=7, bs="ps") + s(Hour, k=24, bs="cc") + s(Hour, Day), data = df_ID1, family = gaussian)
summary(gam3_ID1)$r.sq #0.5616
## [1] 0.5615901
summary(gam3_ID1)$s.table # s(Hour, Day) non significatif
## edf Ref.df F p-value
## s(Week) 3.441641 4.174854 9.9012332 5.390910e-08
## s(Day) 2.134411 2.555628 4.1240830 9.768169e-03
## s(Hour) 21.928036 22.000000 264.7467137 0.000000e+00
## s(Hour,Day) 1.000015 1.000030 0.2120357 6.452041e-01
gam3bis_ID1 <- gam(JOB_DURATION ~ TAP_VERSION + s(Week, k=30, bs="ps") + s(Day, k=7, bs="ps") + s(Hour, k=24, bs="cc") + s(Week, Day), data = df_ID1, family = gaussian)
summary(gam3bis_ID1)$r.sq #0.5617
## [1] 0.5616596
summary(gam3bis_ID1)$s.table # s(Week, Day) non significatif
## edf Ref.df F p-value
## s(Week) 3.441080e+00 4.174153 9.900070e+00 5.408902e-08
## s(Day) 2.134446e+00 2.555668 4.123968e+00 9.766006e-03
## s(Hour) 2.192841e+01 22.000000 2.879624e+02 0.000000e+00
## s(Week,Day) 1.049490e-05 27.000000 2.885171e-07 4.073018e-01
gam3ter_ID1 <- gam(JOB_DURATION ~ TAP_VERSION + s(Week, k=30, bs="ps") + s(Day, k=7, bs="ps") + s(Hour, k=24, bs="cc") + s(Week, Hour), data = df_ID1, family = gaussian)
summary(gam3ter_ID1)$r.sq #0.5651
## [1] 0.5651265
summary(gam3ter_ID1)$s.table # s(Week) non significatif mais s(Week, Hour) significatif
## edf Ref.df F p-value
## s(Week) 1.000001 1.000001 0.2090903 6.475018e-01
## s(Day) 2.118505 2.537928 4.1278985 9.848551e-03
## s(Hour) 21.928301 22.000000 222.5610214 0.000000e+00
## s(Week,Hour) 17.121452 21.307252 2.8290703 1.428118e-05Les combinaisons s(Hour, Day) et s(Week, Day) ne sont pas concluantes car non significatives, mais l’effet bivarié des variables ‘Week’ et ‘Hour’ est lui pertinent même si la fonction de la variable ‘Week’ seule n’est elle pas significative. Nous voyons aussi que la qualité du modèle ne s’est pas améliorée par rapport au modèle 2 puisque nous n’avons pu définir aucun paramètre pour l’effet bivarié. Le meilleur compromis entre ces 3 modèles estimé est donc le dernier pour lequel nous observons ci-dessous son ajustement aux valeurs réelles de Y.
Le constat est le même que pour le modèle 2 ; les faibles volatilités de la durée des JOBS sont mal captées par le modèle. Ainsi, l’ajout de l’effet bivarié n’a pas eu d’impact significatif sur la prédiction comme on peut le constater à travers les \(R^2\) et le graphique ci-dessus qui sont très proches pour les modèles 2 et 3.
4è GAM
Par rapport à une simple fonction s(), les tensors products (te) ont l’avantage de pénaliser les variables dans des directions différentes pour chaque terme, ce qui permet de traiter des variables qui ont des mesures différentes. On peut aussi spécifier les paramètres pour chaque variable (le nombre et le type de splines). De plus, étant donné que la fonction d’interaction inclut automatiquement les fonctions des variables séparément, il n’est pas nécessaire de les rajouter individuellement dans la spécification du modèle.
## [1] 0.8620482
## edf Ref.df F p-value
## s(Day) 3.542003 5.0000 8.626999 3.654141e-10
## te(Week,Hour) 407.030753 443.9237 70.046469 0.000000e+00
Le modèle incluant à la fois les effets bivariés et univarié est la spécification qui explique le mieux le temps des jobs avec une qualité d’ajustement très élevée atteignant 86.2%. De plus, les variables sont toutes significatives ce qui n’était pas le cas des modèles précédents. Le package mgcv propose plusieurs fonctions d’interactions (te, ti et t2) et nous avons fait le choix de tester la fonction ‘t2’ en plus du ‘te’ qui permet de pénaliser avec des fonctions de pénalisation différentes, mais les résultats étaient quasiment les mêmes que pour ce dernier. Sur le plot on voit que les pics sont bien captés par le modèle, comme les fluctuations entre 2 journées et la légère hausse après le pic donc à 5h41 est elle aussi bien modélisée.
Comparaison des modèles
Bien qu’il soit évident que le quatrième modèle généralisé soit le meilleur en termes de qualité prévisionnelle (visuellement et statistiquement: \(R^2\)), nous allons tout de même le vérifier au travers des indicateurs de comparaison présentés précédemment.
| GAM.1 | GAM.2 | GAM.3 | GAM.4 | |
|---|---|---|---|---|
| Nombre d’arguments | 4 | 4 | 5 | 3 |
| Nombre d’arguments significatifs | 3 | 4 | 4 | 3 |
| R2 | 0.3448 | 0.5617 | 0.5651 | 0.862 |
| GCV | / | 21.9 | 21.79 | 7.53 |
| AIC | 31544.84 | 29575.98 | 29550.89 | 24169.54 |
| MAPE | 0.2175 | 0.1483 | 0.1508 | 0.1 |
| Validation | X | ✓ | X | ✓ |
La table ci-dessus nous montre que sur les 4 modèles 2 ont tous les termes significatifs, le modèle n°4 est celui qui minimise les erreurs (GCV, AIC et MAPE) en comparaison à tous les autres ; c’est donc de loin le meilleur modèle pour l’estimation des jobs de la chronique 1. Étant donné que les GAM doivent valider certaines hypothèses sur les résidus, nous allons les vérifier dans la section suivante pour le modèle gardé.
Analyse du meilleur modèle
- C’est la fonction ‘gam.check()’ qui permet la visualisation de la distribution des résidus : la droite de Henry, l’histogramme de distribution et le nuage de points des résidus. On voit alors que nos résidus ne sont pas optimaux mais restent acceptables. Sur la droite de Henry (QQ plot) on voit que les observations sont globalement alignées sur la droite exceptées quelques une en haut à droite du graphique, observations que l’on retrouve en dehors de la ligne à 0 sur le graphique ‘Redids vs. linear pred’ (cf. bulles de la tortue) ainsi que sur le graphique opposant les valeurs prédites aux valeurs observées. On constate aussi 2 groupes principaux dans les résidus, liés au fait que les durées des jobs en journée soient autour de 10/20 secondes puis autour de 40/50 secondes une fois chaque journée.
# Résidus meilleur modèle
par(mfrow=c(2,2))
gam.check(gam4_ID1)
##
## Method: GCV Optimizer: magic
## Smoothing parameter selection converged after 20 iterations.
## The RMS GCV score gradient at convergence was 2.52454e-05 .
## The Hessian was positive definite.
## Model rank = 701 / 701
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## s(Day) 6.00 3.54 1.06 1.00
## te(Week,Hour) 689.00 407.03 0.93 0.04 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
- De même, il est important de vérifier la concurvité des variables, qui correspond à une généralisation de la colinéarité : pour valider l’absence de multicolinéarité les coefficients doivent être compris entre 0 et 1, c’est effectivement le cas pour notre modèle final.
concurvity(gam4_ID1)
## para s(Day) te(Week,Hour)
## worst 0.9780102 0.6595535 5.912765e+28
## observed 0.9780102 0.3857351 2.750745e-03
## estimate 0.9780102 0.3729913 2.039991e-05- Ce plot en 3D permet de voir comment les heures de la journée et les semaines de l’année affectent la durée des requêtes Stitch et donc noter certaines patterns. On retrouve le temps élevé d’éxécution des jobs prédit en début de journée puis quelques fluctuations avec notamment des hausses de la durée autour de 6h ainsi que 15h. Rappelons que les valeurs de la variable ‘Week’ vont de 40 à 52 puis de 1 à 18 donc on retrouve un vide entre ces 2 périodes et quelques fluctuations en fin d’année 2019 et début d’année 2020. On ne trouve pas vraiment de comportement régulier pour la variable ‘Week’. Par ailleurs, les données étant relativement différentes d’une connection ID à l’autre nous devrons réitérer le processus pour chaque chronique.
# Plot 3D des valeurs de Y prédites par le GAM final
# on remodélise en ne prenant que 2 arguments pour faire un graphique 3D
gam_Vis <- gam(JOB_DURATION ~ te(Hour, Week, k=c(24,30), bs = c("cc", "ps")), data = df_ID1, family = gaussian)
# visualisation 3D
vis.gam(gam_Vis, n.grid = 50, theta = 35, phi = 32, zlab = "",
ticktype = "detailed", color = "heat", main = "Plot 3D du modèle GAM N°4")
Autocorrélation
Afin de valider les modèles GAM, les résidus doivent être indépendants et identiquement distribués (iid). Dans le cas des séries temporelles, c’est une hypothèse très forte sur les résidus qui n’est bien évidemment pas respectée en général. Les valeurs des séries temporelles actuelles sont fortement corrélées aux valeurs passées, de sorte que les erreurs du modèle seront également corrélées, c’est ce que l’on appelle le phénomène d’autocorrélation. Cela implique que les coefficients de régression et les résidus estimés d’un modèle peuvent être biaisés. Pour résoudre ce problème, nous pouvons le faire en incluant un processus autorégressif (AR) pour capter l’autocorrélation des erreurs de notre modèle. Nous allons utiliser la fonction gamm() du package mgcv pour réaliser cette nouvelle modélisation du GAM.4 mais cette fois-ci en prenant en compte le phénomène d’autocorrélation.
gam4_ar0 <- gamm(JOB_DURATION ~ TAP_VERSION + s(Day, k=7, bs="cp") + te(Week, Hour, k=c(30,24), bs = c("ps", "cc")),
data = df_ID1,
method = "REML")
gam4_ar1 <- gamm(JOB_DURATION ~ TAP_VERSION + s(Day, k=7, bs="cp") + te(Week, Hour, k=c(30,24), bs = c("ps", "cc")),
data = df_ID1,
correlation = corARMA(form = ~ 1|Week, p = 1),
method = "REML")Nous pouvons comparer les modélisations, sans et avec prise en compte de l’autocorrélation des résidus avec une ANOVA qui sélectionne le meilleur modèle sur la base du log de vraisemblance généralisé. Les résultats du test ci-dessous nous confirment qu’il était judicieux d’estimer en prenant en compte une possible autocorrélation des résidus de la variable ‘Week’. Cela se voit grâce au BIC inferieur pour le gam4_ar1 ainsi que la p-value significative dans la comparaison des 2 modèles.
## Model df AIC BIC logLik Test L.Ratio p-value
## gam4_ar0$lme 1 11 25382.95 25454.60 -12680.47
## gam4_ar1$lme 2 12 25363.18 25441.35 -12669.59 1 vs 2 21.76292 <.0001Afin d’observer l’autocorrélation partielle des résidus normalisés il est possible de réaliser le plot suivant :
Ces graphiques nous donnent les corrélations partielles des résidus avec ses propres valeurs laguées. Les valeurs optimales du pACF se trouvent normalement en dessous des lignes en pointillé ce qui n’est pas vraiment le cas ici, même pour le modèle gam4_ar1. Nous allons donc tenter d’améliorer cela.
Pour choisir les paramètres optimaux du nombre de retards à inclure dans le processus autorégressif sur les résidus, il est possible de passer par la fonction auto.arima du package forecast. Cette fonction choisira le meilleur nombre de retards sur la base de l’AIC. Le code suivant permet de faire cela.
## ar1 ma1 ma2 ma3
## 0.92607777 -0.99713982 0.11062114 -0.02357678D’après la fonction auto.arima, il serait plus pertinent d’appliquer un ARMA(1,3) sur les résidus pour résoudre au mieux l’autocorrélation. Nous devrions vérifier cela en réestimant le modèle et comparant les autocorrélations partielles entre gam4_ar1 et gam4_arma13, mais étant donné le temps conséquent d’exécution, nous avons décidé de présenter les codes sans les exécuter.
#gam4_arma13 <- gamm(JOB_DURATION ~ TAP_VERSION + s(Day, k=7, bs="cp") + te(Week, Hour, k=c(30,24), bs = c("ps", "cc")),
# data = df_ID1,
# correlation = corARMA(form = ~ 1|Week, p = 1, q = 3),
# method = "fREML")
#anova(gam4_ar1$lme, gam4_arma13$lme)
#layout(matrix(1:2, ncol = 2))
#pacf(resid(gam4_ar1$lme, type = "normalized"), lag.max = 48, main = "pACF of gam4 with AR(1)")
#pacf(resid(gam4_arma13$lme, type = "normalized"), lag.max = 48, main = "pACF of gam4 with ARMA(1,3)")ID 2
1er GAM
Nous démarrons les modélisations pour la chronique n°2 avec un modèle GAM basique où nous laissons le modèle ajuster ses paramètres seul, c’est à dire avec la méthode “REML”. Pour rappel les connections ID égales à 2 sont observées l’été 2019, la variable ‘TAP_VERSION’ prend seulement 2 valeurs (‘1.4.33’ pour 2366 observations et ‘1.4.34’ pour les restantes). Pour cette sous-base de la chronique n°2 la valeur des heures considère aussi les minutes que nous avons codées comme étant les décimales après la virgule pour pouvoir passer la variable en quantitative précédemment. Dans un premier modèle généralisé nous incluons donc la variable ‘TAP_VERSION’ bien qu’elle ne prenne que 2 valeurs nous supposons qu’elle ne sera pas significative mais cherchons à en être sûrs pour ne pas passer à côté d’une éventuelle information. Nous ajoutons aussi des fonctions des variables ‘Week’, ‘Day’ et ‘Hour’.
##
## Family: gaussian
## Link function: identity
##
## Formula:
## JOB_DURATION ~ TAP_VERSION + s(Week, k = 1) + s(Day, k = 1) +
## s(Hour)
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 178.204 0.419 425.309 <2e-16 ***
## TAP_VERSION1.4.34 -1.626 2.014 -0.807 0.42
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(Week) 1.010 1.020 36.42 1.53e-09 ***
## s(Day) 1.951 1.998 20.89 9.87e-10 ***
## s(Hour) 8.721 8.977 45.47 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.164 Deviance explained = 16.8%
## -REML = 11072 Scale est. = 410.5 n = 2497
On voit d’emblée que la variable ‘TAP_VERSION’ n’est pas significative comme nous l’avons supposé précédemment, mais toutes les fonctions incluses elles le sont. De plus, on remarque que l’évolution des variables de semaine et de jours est inverse par rapport aux estimations de l’ID1. La variable ‘Week’ qui était incurvée est désormais parfaitement linéaire (mais on observe tout de même que les semaines fin août ont un temps d’éxécution plus long que les semaines début juillet avec une hausse linéaire) et celle des jours ‘Day’ est maintenant légèrement concave. Leurs degrés de liberté estimés confirment cela car ils sont proches de 1, nous garderons donc pour les estimations suivantes que les variables ‘Day’ et ‘Hour’. Concernant cette dernière on voit que les temps des jobs ne suivent pas les mêmes évolutions dans une journée que l’ID1, ce que nous avions pu noter dès la visualisation des données dans une partie antérieure, et l’on voit que le pic est désormais autour de 9-10h du matin et que l’après-midi les temps augmentent à nouveau pour diminuer en fin de journée. Finalement, on voit que la qualité d’ajustement du modèle est vraiment faible (-20 points par rapport aux GAM de la chronique une) ; le \(R^2\) est de 16.4% pour ce premier modèle, nous cherchons à augmenter le plus possible cette qualité et résumerons notre recherche en un tableau récapitulatif.
Visuellement la précision faible de nos prédictions se confirme avec une courbe des valeurs prédites par le modèle qui ne suit pas les évolutions des valeurs réelles et capte très mal le pic de la matinée à 9h09.
Optimisation du modèle
|
Effets univariés seulement
|
Effets bivariés
|
||||||
|---|---|---|---|---|---|---|---|
| GAM.1 | GAM.2 | GAM.3.1 | GAM.3.2 | GAM.3.3 | GAM.3.4 | GAM.3.5 | |
| Termes inclus dans le modèle | |||||||
| 1er terme | TAP_VERSION | s(Day,k=7) | te(Day,Hour,k=c(7,48)) | te(Week,Hour,k=c(7,48)) | te(Day,Hour,k=c(7,48) | te(Day,Hour,k=c(7,48)) | te(Day,Hour,k=c(7,48)) |
| 2è terme | s(Week,k=1) | s(Hour,k=48) | te(Week,Hour,k=c(7,48)) | te(Week,Hour,k=c(7,48)) | |||
| 3è terme | s(Day,k=1) | te(Week,Day,k=c(7,7)) | |||||
| 4è terme | s(Hour) | ||||||
| Nombre de termes significatifs | 3/4 | 2/2 | 1/1 | 1/1 | 1/1 | 2/2 | 3/3 |
| R2 | 0.164 | 0.3576 | 0.3779 | 0.3742 | 0.3722 | 0.4204 | 0.4313 |
| GCV | / | 321.82 | 317.5 | 317.82 | 319.8 | 302.04 | 296.87 |
| AIC | 22127.24 | 21504.82 | 21468.33 | 21471.76 | 21486.68 | 21338.68 | 21295.04 |
| MAPE | 0.0622 | 0.0504 | 0.0488 | 0.0502 | 0.049 | 0.0459 | 0.0457 |
| Validation | X | ✓ | ✓ | ✓ | ✓ | ✓ | ✓ |
Le tableau ci-dessus résume les différentes modélisations que nous avons mises en place pour trouver les meilleures prédictions du temps des jobs sur la chronique 2. Le modèle n°2 reprend les variables significatives du premier GAM en spécifiant le nombre de splines appliquées; 7 pour les jours de la semaine et 48 pour les heures car comme précisé avant, chaque heure comporte 2 observations. La qualité grimpe de 16.4% à 35.76% soit presque 20 points de plus, puis dans un troisième GAM nous ajoutons les effets bivariés que nous appliquons pour chaque variable. Il apparaît alors que les tensors products permettent d’améliorer grandement le modèle et le \(R^2\) le plus élevé que nous obtenons est de 43.13% correspondant au modèle “3.5” où nous avons inclu les effets bivariés des jours/heures, semaines/heures et semaines/jours. En plus de maximiser le \(R^2\), ce modèle minimise l’erreur par CV (GCV), le critère Akaike et le MAPE qui est de 4.57% d’erreur absolue moyenne, ce qui est plus faible que pour les modélisations de l’ID 1 qui variaient entre 10 et 15%. Nous analysons à présent le modèle 3.5 sélectionné.
Analyse meilleur modèle
- Commençons par confronter visuellement les valeurs prédites aux valeurs réelles de “JOB_DURATION” du modèle 3.5 gardé. Nous voyons alors sur le graphique ci-dessous que les variations sont beaucoup mieux captées que le modèle initial où la méthode ‘REML’ était spécifiée dans le but que les paramètres soient éstimés seuls. On voit combien il est important de spécifier le nombre de splines et de considérer des effets bivariés entre les variables. Ainsi, sur une semaine prise au hasard dans notre échantillon d’apprentissage, on voit que les valeurs prédites par le modèle sont assez proches de la réalité.
- De la même manière que pour le modèle retenu pour la chronique une nous regardons la distribution des résidus. Nous observons alors comme visible ci-après, que certains résidus se distinguent des autres puisqu’ils correspondent à des temps de jobs importants (entre 300 et 400 secondes d’éxecution) mais prédits à la baisse par le modèle (entre 180 et 220 secondes). Ces quelques points sont facilement identifiables sur les 4 graphiques de diagnostic des résidus.
##
## Method: GCV Optimizer: magic
## Smoothing parameter selection converged after 13 iterations.
## The RMS GCV score gradient at convergence was 0.0008146698 .
## The Hessian was positive definite.
## Model rank = 654 / 654
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## te(Day,Hour) 335.0 75.6 1.06 0.90
## te(Week,Hour) 276.0 44.4 0.99 0.42
## te(Week,Day) 42.0 28.2 0.98 0.29
- Concernant la concurvité on constate là encore qu’il n’y a pas de problème de multicolinéarité car les coefficients de concurvité observés prenent des valeurs entre 0.4 et 0.85 pour les 3 termes du modèle final.
## para te(Day,Hour) te(Week,Hour) te(Week,Day)
## worst 1.631445e-12 1.0000000 1.0000000 1.0000000
## observed 1.631445e-12 0.8466615 0.4007018 0.5302344
## estimate 1.631445e-12 0.1565188 0.1952045 0.2851099- Enfin, nous voyons l’évolution des durées des requêtes prédites par le GAM retenu, selon les variables ‘Day’ et ‘Hour’. Nous regardons leur effet car par rapport aux autres variables avec effet bivarié, ce sont celles qui ont le plus grand effet sur Y (cf. summary du modèle). Les temps des jobs sont plus importants en milieu de semaine, on retrouve toujours le pic de début de journée ainsi que celui en fin d’après-midi. Le pic important est une pattern commune aux temps des jobs de la chronique 1, sauf que l’un est à 5h41 tous les jours tandis que l’autre est à 9h09 - ce qui peut être lié à la période étudiée entre l’année scolaire et les périodes de grandes vacances.
Autocorrélation
De la même façon que pour l’ID1, nous pouvons observer la pertinence de rajouter à l’estimation des résidus du modèle GAM un processus autorégressif d’ordre 1. Le résultat du test Anova ci-dessous confirme encore une fois la nécéssité de prendre en compte le phénomène d’autocorrélation dans l’estimation de l’ID2 par le modèle GAM 3.5. Cela s’observe à travers un AIC et un BIC plus faibles pour le modèle gam3.5_ar1 et par le fait que cette différence soit significative car la p-value est très proche de 0.
## Model df AIC BIC logLik Test L.Ratio p-value
## gam3.5_ar0$lme 1 11 21555.53 21619.57 -10766.77
## gam3.5_ar1$lme 2 12 21484.09 21553.94 -10730.04 1 vs 2 73.44542 <.0001Les plots ci-dessous nous montrent l’autocorélation partielle des résidus normalisés:
La prise en compte de l’autocorrélation n’est pas très convaincante dans le graphique de droite. Il faudrait que les traits verticaux ne dépassent pas le seuil en poitillé bleu pour que les résidus soient considérés comme non-autocorrélés. La solution serait alors d’optimiser les paramètres du modèle autorégressif en augmentant l’ordre du comportement AR ou en ajoutant une partie MA au processus. Dans notre analyse, il nous a semblé nécéssaire d’expliquer comment traiter l’autocorrélation, mais n’avons pas optimisé la prise en compte de l’autoccorélation pour chaque chronique car le principe est le même et le temps d’exécution est très long.
ID 3
Optimisation
Enfin, c’est avec la chronique correspondant à la connection ID 3 que nous terminons l’analyse et la construction de modèles GAM. Comme observé en partie de visualisation des données, en plus des tendances au sein d’une journée nous observons pour cette chronique des tendances sur une semaine puisqu’on ne retrouve pas de pics de temps plus importants dans la journée le samedi. Nous incluons donc toutes les variables et construisons les modèles en les améliorant au fur et à mesure les prédictions, comme pour les estimations des chroniques 1 et 2.
|
Effets univariés seulement
|
Effets bivariés
|
||
|---|---|---|---|
| GAM.1 | GAM.2 | GAM.3.1 | |
| Termes inclus dans le modèle | |||
| 1er terme | TAP_VERSION | TAP_EXIT_STATUS | TAP_EXIT_STATUS |
| 2è terme | TAP_EXIT_STATUS | s(Week,k=10) | te(Week,Day,k=c(10,7)) |
| 3è terme | s(Week) | s(Day,k=7) | te(Hour,Day,k=c(24,7)) |
| 4è terme | s(Day,k=1) | s(Hour,k=24) | te(Hour,Week,k=c(24,10)) |
| 5è terme | s(Hour) | ||
| Nombre de termes significatifs | 4/5 | 4/4 | 4/4 |
| R2 | 0.639 | 0.681 | 0.885 |
| GCV |
|
62305.42 | 26049 |
| AIC | 22149.91 | 21969.55 | 20544.62 |
| MAPE | 0.104 | 0.0969 | 0.0445 |
| Validation | X | ✓ | ✓ |
Comme visible sur la table ci-dessus nous avons estimé 3 modèles pour cette dernière sous-base de notre échantillon d’entraînement. Dans le premier nous incluons toutes les variables potentielles pour expliquer le temps des jobs de CONNCETION_ID=3 : il en ressort que ‘TAP_VERSION’ n’est pas significative donc nous l’excluons des estimations. Nous voyons aussi que contrairement aux ID 1 et 2, la qualité d’ajustement est directement élevée (63.9%) donc les variations doivent être plus régulières pour que le modèle les capte mieux. Nous améliorons cette première estimation dans une deuxième GAM où nous spécifions le nombre de splines ; de la même manière que les estimations précédentes nous fixons k au maximum c’est-à-dire au nombre de valeurs que prend chaque variable (dans la 3è sous-base il y a 10 semaines). Aussi, comme les modélisations précédentes, le fait de spécifier les types de splines n’améliore pas les estimations voire les diminue (très légèrement). En outre, dans un dernier modèle nous ajoutons les effets bivariés avec les 3 combinaisons possibles des variables explicatives et atteigons alors un \(R^2\) de 0.885, soit le meilleur que nous ayons pu trouver jusqu’ici. Regardons plus précisément les caractéristiques de ce modèle n°3 que nous gardons comme modèle final.
Analyse meilleur modèle
- Visuellement nous constatons la bonne qualité des prévisions réalisées à partir du dernier modèle ; les valeurs prédites suivent les évolutions principales des valeurs observées, la seule différence résidant dans quelques temps des jobs qui sont prédits pour la plupart plus faibles que la vraie valeur.
- Le diagnostic montre cette fois des résidus beaucoup mieux distribués que pour les modèles des ID 1 et 2. Leur distribution visible sur l’histogramme s’apparente davantage à celle d’une loi normale, sur le graphique “Resids vs. linear pred.” les valeurs sont plus groupées autour de 0 et sur le graphique des valeurs prédites et observées elles sont davantage situées sur la bissectrice. Les résidus de ce troisième modèle sont donc bien meilleurs que les précédents.
##
## Method: GCV Optimizer: magic
## Smoothing parameter selection converged after 18 iterations.
## The RMS GCV score gradient at convergence was 0.001604292 .
## The Hessian was positive definite.
## Model rank = 439 / 439
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## te(Week,Day) 69.0 62.1 1.00 0.54
## te(Hour,Day) 161.0 134.6 1.04 0.87
## te(Hour,Week) 207.0 54.1 1.00 0.50
- Il n’y a pas non plus de problème de colinéarité entre les termes de l’estimation puisque les coefficients de concurvité se situent pour ce modèle entre 0.7 et 0.98.
## para te(Week,Day) te(Hour,Day) te(Hour,Week)
## worst 0.3631692 1.841721e+28 0.9993514 0.9987963
## observed 0.3631692 9.773962e-01 0.8631211 0.7027251
## estimate 0.3631692 2.204587e-01 0.1704511 0.1422472- Finalement, nous pouvons voir le plot 3D des valeurs de Y prédites par le modèle additif généralisé final ; nous regardons les variations des temps des jobs en mettant 2 à 2 les variables pour lesquelles nous avons supposé un effet bivarié. À partir du premier plot on voit que le temps des jobs est important surtout en fin d’après-midi et en milieu de soirée, et est plus faible pour les premières semaines (début janvier) et les dernières (mi-mars) comparé aux semaines entre ces deux périodes. Par la variable ‘Day’ nous constatons que les temps d’éxécution des requêtes sont plus courts en fin de semaine, c’est à dire le week end et dans la nuit du dimanche à lundi.
Pour conclure nous pouvons voir à la fois par le diagnostique des réisuds, les différents graphiques et les taux d’erreur que le modèle estimé pour cette troisième sous-base est le plus précis que nous ayons pu construire jusqu’ici.
Autocorrélation
Encore un fois nous pouvons observer les critères d’information de l’AIC et du BIC plus faibles significativement pour le modèle gam3_ar1, c’est à dire où est prise en compte l’autocorrélation des résidus.
## Model df AIC BIC logLik Test L.Ratio p-value
## gam3_ar0$lme 1 12 20989.41 21053.78 -10482.71
## gam3_ar1$lme 2 13 20726.02 20795.75 -10350.01 1 vs 2 265.3959 <.0001Les plots ci-dessous nous montrent l’autocorélation partielle des résidus normalisés :
Pour la chronique 3, nous pouvons observer visuellement une amélioration des résidus lors de la prise en compte de l’autocorrélation par un AR d’ordre 1. En effet les autocorrélations partielles des résidus normalisés sont majoritairement comprises entre les lignes en pointillé. Malgré cela, il serait tout de même judicieux d’optimiser l’ordre du processus AR pour une meilleure prise en compte de ce phénomène comme nous l’avons développé dans la chronique 1.
2.3 Résumé visuel des modèles retenus
Pour terminer cette partie de recherche et d’optimisation des modèles additifs généralisés, nous regardons les graphiques d’évolution des valeurs observées et prédites par chaque modèle retenu.
3-Soumission des prévisions et conclusion
Maintenant que nous avons construit les modèles les plus précis en termes de prévisions nous allons soumettre des prédictions pour la compétition kaggle. Pour cela nous devons mettre la base ‘test’ au même format que la base ‘train’ pour pouvoir appliquer les modèles construits à partir de cette dernière. La base test contient 1700 observations et 5 variables comme précisé en introduction : notre objectif est donc, à partir des données disponibles, de prédire le temps d’éxécution des jobs Stitch.
En reprenant les commandes pour la base d’apprentissage, nous séparons la base test en 3 sous-échantillons à partir de la CONNECTION ID, puis nous créons les variables de semaine, jours et heures pour chacune d’entre elles grâce à la variable ‘START_TIME’. Enfin, nous passons les variables créées au format numérique pour pouvoir appliquer librement les modèles additifs généralisés.
Nous commençons par appliquer le quatrième modèle sélectionné pour la sous-base de la chronique 1, grâce à la fonction predict.gam() contenue elle aussi dans le package mgcv sous R. Nous réalisons ainsi les prévisions pour les 3 base de données test en créant chaque fois une colonne ‘Predicted’ dans l’échantillon concernée, puis nous supprimerons toutes les variables autres que celle-ci et l’ID, pour retrouver la forme du fichier sample_solution que nous importerons sur Kaggle pour voir les résultats de notre estimation, après avoir rassemblé de nouveau les 3 sous-bases.
Les lignes de code suivantes ont permis un tel processus :
#------- ID 1 ----------
# On fait les prévisions pour la base test de la chronique 1 avec le modèle sélectionné pour cet ID
#test_ID1$Predicted <- predict.gam(gam4_ID1, test_ID1) #ne fonctionne pas car les versions TAP sont différentes, nous retirons cette variable et commencons une nouvelle fois
# On reestime le modèle 4 en enlevant 'TAP_VERSION'
gam4_ID1 <- gam(JOB_DURATION ~ s(Day, k=7, bs="cp") + te(Week, Hour, k=c(30,24), bs = c("ps", "cc")), data = df_ID1)
# on fait les prévisions
test_ID1$Predicted <- predict.gam(gam4_ID1, test_ID1)
#------- ID 2 ----------
test_ID2$Predicted <- predict.gam(gam3.5_ID2, test_ID2)
#------- ID 3 ----------
test_ID3$Predicted <- predict.gam(gam3_ID3, test_ID3)
# On remet les bases au format du fichier 'sample_solution'
# on retire toutes les autres variables
test_ID1 <- test_ID1[,-c(1:4,6:8)]
test_ID2 <- test_ID2[,-c(1:4,6:8)]
test_ID3 <- test_ID3[,-c(1:4,6:8)]
# on regroupe les 3 bases un 1 df à soumettre sur kaggle
sample_solution <- rbind(test_ID1, test_ID2, test_ID3)
# on exporte ces 2 df au format csv
rio::export(sample_solution, "./sample_solution.csv")Après soumission, nous avons obtenons un taux d’erreur MAPE de 14.28%. Notre étude est donc terminée et les prévisions du temps d’éxecution des requêtes Stitch sont satisfaisantes.